d-wave pairing in lightly doped Mott insulators 
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We define a suitable quantity Z c that measures the pairing strength of two electrons added to the 
ground state wave function by means of the anomalous part of the one-particle Green's function. 
Z c discriminates between systems described by one-electron states, like ordinary metals and band 
insulators, for which Z c = 0, and systems where the single particle picture does not hold, like 
superconductors and resonating valence bond insulators, for which Z c 7^ 0. By using a numerically 
exact projection technique for the Hubbard model at U/t = 4, a finite value of Z c , with d-wave 
symmetrry, is found at half filling and in the lightly doped regime, thus emphasizing a qualitatively 
new feature coming from electronic correlation. 

PACS numbers: 71.10.-w, 71.10.Fd, 71.27.+a 
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Since the discovery of high-temperature (HTc) super- 
conductors, the question whether a strongly correlated 
system, containing only repulsive interactions, may dis- 
play a superconducting ground state, has been inten- 
sively debated until now, mainly because HTc super- 
conductors arc certainly strongly correlated and their 
critical temperature cannot be explained by a standard 
electron-phonon mechanism. Indeed, any mean-field ap- 
proximation fails to explain superconductivity without 
an electron-electron attraction, that could be either ex- 
plicit, like in the attractive Hubbard model, or mediated 
by some boson, like in the BCS theory. On the other 
hand, it is well known that correlated wave functions, 
that highly improve the mean-field energies, usually favor 
superconductivity. 01 S S S S] Recently, the approach 
based on projected superconducting wave functions has 
been recently renewed, as it is possible to reproduce many 
important experimental aspects just by assuming that 
the ground state is described by a projected BCS wave 
function. 0, S| Within this scheme, the superconductiv- 
ity is "hidden" in the insulating state, where phase co- 
herence is inhibited by the charge gap, and it is indeed 
stabilized by a small amount of doping. Unfortunately, 
the validity of this scenario for the actual ground state 
of a microscopic model remains an highly debated and 
controversial issue, despite a huge amount of numerical 

work. 10, Si, Eg 

In this letter, we present a numerical study in favor 
of the mentioned scenario. We consider the single-band 
Hubbard model on L-site clusters: 



H=-t £ 



h.c. + Uy]nj t fiii,i - fJ-N, (1) 



where c\ a creates an electron with spin a at the site 
i, ni,<j = c\ a Ci,a is the density operator at the site i, 
fx is the chemical potential, and TV = a ^ s ^ ne 
total number of particles. We use square lattices tilted 
by 45° with L 2 = 2l 2 and / odd with periodic boundary 
conditions, so that the half-filled case (N = L) has a non- 
degenerate ground state even at U = 0. This condition 



strongly reduces the size effects, particularly important 
in two dimensions (2D). 

Our purpose is to study the anomalous part of the 
equal-time Green's function at zero temperature: 
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where \^o) is the ground state of the Hamiltonian H. 
Obviously, F k BSP can be non zero only in the thermo- 
dynamic limit, where, because of the broken symmetry 
phenomenon (BSP), |^o) may not have a definite num- 
ber of particles. On a finite system BSP does not occur 
and F^ SP is always zero, the ground state V^o) having 
always a definite number of particles N. Therefore, on 
finite systems, it is useful to consider a related quantity: 
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where ) is the ground state with TV particles. Fol- 
lowing well known concepts of BSP, by adding to the 
Hamiltonian a small anomalous term, that breaks the 
symmetry, proportional to 5, and taking the limit 5 — ► 0, 
after the thermodynamic limit, it is easy to convince that 
the finite-size expression F k converges to F PSP . 

In order to evaluate Fk, we consider a trial BCS wave 
function |^g) containing only components with even 
number of particles: 
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where |0) is the vacuum, f k = A k /(e k - /Lt + E k ), 
e k = —2t(cosk x + cos k y ) is the free-electron dispersion, 
Ho is a variational parameter playing the role of the chem- 
ical potential of the BCS Hamiltonian, A k is the corre- 
sponding gap function chosen to have either a d-wave 
symmetry A k = A(cosfca; — cosk y ) or an s-wave sym- 
metry A fe = A, and E k = ^(e k - fi ) 2 + A' 2 . The cor- 
relation factor J = e~ 9 *-'i ni ^ ni > 1 , g being a variational 
parameter, partially projects out expensive energy con- 
figurations with doubly occupied sites. In the following, 
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we assume that, as verified on small clusters, for an ap- 
propriate gap function symmetry, the TV-particle ground 
state of H have non-vanishing overlap with |^g)> a t least 
for N = N* and N = N* + 2, where N = N* represents 
a closed shell density. Then, on each finite system, we 
compute for fixed imaginary time r: 



(* G | e -^ /2 iAa e ~ Tg/2 l^) 



(5) 



Fk(i~) will be finite even on a finite size, because the trial 
function |^g) contains sectors with different N, namely 
I*g) = E w inI^g)- Then for lar g° r: 



F k (r) = F k 



~>N*+2 c -tA c /2 



b N * 



(6) 



where N* is chosen by fixing the true chemical potential 
/i, A c is the energy gap between the states with N* and 
N* + 2, and b N = (^^\^^)a N . Therefore, in order to 
estimate F k , the overall factor 7 = bN* + 2/bN*e~ T ^ a ^ 2 
of Eq. 10 must be computed by evaluating the average 
number of particles (N) = A* + 27 2 /(l+7 2 ). An efficient 
choice of \^g) is obtained by tuning /^o very close to the 
energy level e k with N* + 2 particles, and A <C \e k — Mo|- 
In this limit, \^g) essentially contains only two compo- 
nents in the sectors with N = N* and N = N* + 2 and 
the imaginary time projection can be done without pay- 
ing much attention to the chemical potential, the number 
of particles being conserved by H . 

The superconducting order parameter P is simply 
related to the short-distance component of the real- 
space Fourier transform Fr = l/i^e'^Ft, namely 
P d = 2F V , with Tj = (±1,0) or 77 = (0, ±1) for a d- 
wave superconductor with A k = A(cosfc x — cosfcj,), and 
P s = F n with r; = (0, 0) for an s-wave superconductor 
with Afe = A. Based on the resonating-valence bond 
(RVB) theory, we want to define a quantity, closely 
related to F k , that makes sense also when Pd = 0, as in 
the insulating case. Indeed, according to the RVB theory, 
the insulator already contains some sort of pairing, the 
electrons being paired in RVB singlets. Following this 
paradigm, we normalize the anomalous pairing function 
in real space and define the quantity: 



Fx 



9r 



(7) 



where the sum is for all distances in the lattice (compat- 
ible with the symmetry of Fr). Notice that, also for an 
infinite system with a charge gap, like the half-filled Hub- 
bard model, though the anomalous average Fr is zero, 
a finite ratio gu is still possible. g\ may be interpreted 
as the probability function of two electrons added into a 
singlet state at a distance R and, therefore, determines 
the pair function, which is the Cooper pair if P > 0. For 
a BCS superconductor with a full gap (e.g., s-wave) gu 
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FIG. 1: Size scaling of the d-wave pairing strength Z c for 
different wave functions at half-filling: the uncorrelated BCS 
wave function with A = 0.05, the free-electron wave func- 
tion (A — > 0), the wave function with a finite antiferromag- 
netic order parameter [obtained by adding to the A — » 
BCS Hamiltonian a further antiferromagnetic mean field term 
(— I) 1 AAF(ni,i — n i,i)]> an d the Gutzwiller wave function with 
g = 1 (full squares). The same quantity for the exact ground 
state of the two-chain t— J ladder at low doping is also shown 
(full triangles). 



is localized with a well defined coherence length, whereas 
if the pairing function has nodes, like in the d-wave case, 
gu decays with a power law, but, nonetheless, it exists 
a finite amplitude of the pairing function at short dis- 
tance. In strong-coupling superconductivity, the coher- 
ence length is expected to be small, and the most impor- 
tant contribution to the pairing function is at the shortest 
distance ?y allowed by the symmetry of the pairing func- 
tion. Therefore, we define a fundamental quantity, the 
pairing strength, that we denote by: 



z c =\g v \- 



(8) 



As shown in Fig. ^ m the thermodynamic limit, Z c 
is finite in an uncorrelated d-wave superconductor and, 
instead, it vanishes for free electrons, for a weakly- 
correlated Fermi liquid, described by a Gutzwiller wave 
function, and for a spin-density wave state, with a finite 
antiferromagnetic order parameter. By introducing the 
pairing strength Z Cl we provide a general and quantita- 
tive definition of pairing, which holds not only for sim- 
ple superconductor systems, but also for non-BCS sys- 
tems. A simple example is the two- leg ladder (see Fig.^l, 
where off-diagonal long-range order is suppressed by one- 
dimensional quantum fluctuations, but, nevertheless, Z c 
remains finite, showing that pairin g is well defined also 
in this system, as widely accepted. |l2l fl3l IT^j 

There is also another important reason to study the 
pairing strength Z c instead of the order parameter P: in 
a strongly correlated system the value of the quasipar- 
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FIG. 2: The d-wave pairing strength Z c (circles) as a function 
of U in an L = 18 site cluster with N* = L. The quasiparticle 
weight Zn for N = 16 (squares) at the momentum (27r/3,0) 
and the short-distance anomalous average Pd (triangles) are 
also shown. Continuous {Zn), long dashed (Z c ) and dotted 
(Pd) lines correspond to the variational calculation with the 
variational wave function of Eq. Q with an optimized 
and projected onto the subspaces with 16 and 18 electrons. 
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ticle weight Z n, d enned as Z 

be very small [la, [l6| and P, if finite, is expected to be 
at most of the same order. Therefore, whenever Zn is 
very small, it is very difficult to detect a non-zero value 
of the anomalous average P. The suppression of d-wave 
pairing obtained in previous calculations 0, Q can be 
explained by the fact that the quasiparticle weight de- 
creases very rapidly with U (see below and Ref. fl7j). 
The question of a finite anomalous average of order Zn 
could be still compatible with the published numerical 
calculations and it is beyond the scope of this work. For 
this reason it is very difficult, at the present time, to de- 
tect pairing by studying directly the order parameter or 
the pairing correlations. On the other hand, the pairing 
strength Z c , being a ratio of two quantities of the same 
order, is not affected by the small quasiparticle weight 
Zjq and, therefore, represents a much more sensitive de- 
tector of superconductivity. Strictly speaking, the fact 
of having a finite Z c in the thermodynamic limit does 
not imply that also P is finite, but only the presence of 
paired electrons, as in the example of the two-leg ladder. 
However, in 2D, whenever the compressibility is finite, 
it is reasonable to expect condensation of pairs. Unfor- 
tunately, this assumption cannot be supported by direct 
calculation in the case of the Hubbard model. Neverthe- 
less, we can safely discuss the d-wave pairing properties 
of the Hubbard model because a non-zero value of Z c 
appears very clearly. 

Before considering the quantum Monte Carlo results, 
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FIG. 3: The d-wave pairing strength Z c as a function of 
the projection time r for different clusters at half-filling. The 
inset shows the finite size scaling of Z c — Z^ for r = 6. 



we show the results for a small lattice of L = 18 sites, 
where it is possible to perform exact diagonalization. In 
Fig. |21 we report the results for Z c , P^ and Zn as a 
function of U/t. The comparison between P^ and Zn 
clearly supports our picture: P^ is small in the strong- 
coupling Hubbard model because Zn is small. On the 
other hand, a much higher signal of paired electrons can 
be obtained by studying the pairing strength Z c , which 
shows a broad maximum Z c ~ 0.5 for U/t ~ 16, de- 
creasing to the minimum value Z c ~ 0.16 for U — > oo. 
Although the value of the anomalous average Pd and the 
quasiparticle weight are highly overestimated at strong 
coupling, the variational wave function of Eq. Q with 
an optimized provides an excellent estimate of Z c for 
all values of U/t < 20, thus capturing the correct feature 
of pairing. 

We now turn to larger systems and calculate the pair- 
ing strength by using the zero-temperature Monte Carlo 
projection technique based on auxiliary fields. [l8j The 
inclusion of the simple Gutzwillcr factor J , which im- 
proves substantially the convergence in imaginary time 
r, is particularly important because, at finite doping, the 
sign problem prevents us to work with arbitrary large 
imaginary time. 

Firstly, we consider the half-filled case. Remarkably, 
for large systems, Z c increases with the projection time r, 
see Fig- El In this case, in order to reduce size effects, wc 
have performed the finite-size scaling of Z c — Zf? , where 
is the value corresponding to the Gutzwiller wave 
function |* G ) with A — * (Zf — > in the thermody- 
namic limit, see Fig. QJ. Already for L = 50 sites, the 
evaluation of the pairing strength is rather accurate and 
close to larger sizes. As shown in Fig. 01 Z c increases 
monotonically with r, and the value for r = 6 should 
safely represent, for all sizes considered, a rather accurate 
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FIG. 4: The d-wave pairing strength Z c — Z^ as a function 
of the projection time r for different dopings. The arrows 
indicate the value of Z c corresponding to the variational wave 
function of Eq. (|1J with an optimized for and 8 holes. 
Z c for the s-wave pairing in one dimension is also shown. 

lower bound for the r-too limit. The thermodynamic 
limit of Z c appears therefore clearly finite, considering 
also that for r = 6 Z c — Z^f increases with the system 
size, see the inset of Fig. |3] Our extrapolation is consis- 
tent with a finite Z c > 0.1 in the thermodynamic limit. It 
is worth noting that, in this finite value of Z c does 

not mean that the ground state is superconducting, but 
only that two added electrons in the half-filled Hubbard 
model are paired together in a d-wave singlet. 

In the doped case the limitation of the sign problem 
and the strong dependence upon hole doping prevent us 
to perform an accurate finite size scaling, but, as shown 
in Fig. the effect remains rather clear, and Z c — 
increases immediately as soon as the projection time is 
turned on. On the other hand, at large enough doping, 
we have found the opposite effect, clearly against d-wave 
pairing. Interestingly, also the s-wave pairing strength, 
obtained by using an initial wave function \^g) with s- 
wave symmetry, is not enhanced by the imaginary time 
projection (not shown). However, we have not checked 
other symmetry sectors like p or d xy and, therefore, at 
large doping there maybe other types of pairing insta- 
bilities. [19j It is remarkable that also for 50 sites the 
simple Gutzwillcr wave function with a finite d-wave gap 
function (stabilized in the undcrdoped region) describes 
accurately the enhanced pairing strength Z c , as shown 
by the arrows in Fig. 

Although we worked at zero temperature, one should 
be tempted to associate Z c with the pseudogap observed 
in the HTc superconductors. In this scenario, the onset 
of the pseudogap, for a temperature T*, is marked by 
a finite value of Z c , indicating paired electrons without 
phase coherence, then at a much lower temperature T c 
the pairs eventually condense, giving rise to a true long- 



range order. [2p| 

In conclusion, it turns out that the pairing is a robust 
property of lightly-doped Mott insulators and appears 
already in small size calculations. The pairing strength 
Z c increases with decreasing doping and has its maxi- 
mum at half filling, where phase coherence is inhibited 
by the charge gap. Notice that this is a peculiar feature 
of the 2D model, as in one dimension no evidence of a 
finite pairing strength is found even at half filling, sec 
Fig. 0] We argue that a finite pairing strength for an 
insulator is just the qualitative new feature that discrim- 
inates a band insulator from a 2D RVB insulator, which 
is defined in terms of singlet pairs, as in a fully projected 
BCS wavefunction of Eq. Q . Within this definition, de- 
termined by the measurable quantity Z Cl and indepen- 
dent of the variational ansatz, we discover the possibility 
to have an RVB-like insulator also when the existence of 
a finite pairing strength is accompanied by the antiferro- 
magnetic long-range order, like in the half-filled Hubbard 
model. [2lJ In this simple model, the insulator is somehow 
prepared to become superconductor (with d-wave sym- 
metry) and a small amount of doping allows the propaga- 
tion of the RVB pairs. This scenario offers a simple 
and natural explanation of the ultimate mechanism of 
HTc superconductivity. 
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